1 Object

  • Summarize what I did for data downloading

2 gridMET

  • Target variables:
    • “pr”: precipitation (mm)
    • “tmmn”: Minimum Near-Surface Air Temperature (Kelvin)
    • “tmmx”: Maximum Near-Surface Air Temperature (Kelvin)
    • “pet” : Reference grass evaportranspiration (mm)
  • Using well-level data (point data), extract values from gridMET raster data (2007-2015)
    • extract only in season daily data (i.e., from 4/1 to 9/30)
  • change those units
    • pr: mm to inch
    • tmmn: K to C
    • tmmx: K to C
    • pet: mm to inch
  • calculate daily GDD from tmmn and tmmx
  • summarize those variables by year and well id

3 SSURGO

3.1 Description about target variables

  • Target variables (reference website: here[https://sdmdataaccess.sc.egov.usda.gov/documents/TablesAndColumnsReport.pdf])
    • sandtotal_r:
      • %
      • Mineral particles 0.05mm to 2.0mm in equivalent diameter as a weight percentage of the less than 2 mm fraction
    • silttotal_r:
      • %
      • Mineral particles 0.002 to 0.05mm in equivalent diameter as a weight percentage of the less than 2.0mm fraction
    • claytotal_r:
      • %
      • Mineral particles less than 0.002mm in equivalent diameter as a weight percentage of the less than 2.0mm fraction
    • ksat_r:
      • um/s
      • The amount of water that would move vertically through a unit area of saturated soil in unit time under unit hydraulic gradient.
    • awc_r:
      • cm/cm
      • The amount of water that an increment of soil depth, inclusive of fragments, can store that is available to plants. AWC is expressed as a volume fraction, and is commonly estimated as the difference between the water contents at 1/10 or 1/3 bar (field capacity) and 15 bars (permanent wilting point) tension and adjusted for salinity, and fragments.
    • slope_r
      • %
      • The difference in elevation between two points, expressed as a percentage of the distance between those points. (SSM)
  • There are many similar variables depends on the particle size
    • Ex) in the sand case
      • sandtotal_ : Mineral particles 0.05mm to 2.0mm in equivalent diameter
      • sandvc_ : Mineral particles 1.0mm to 2.0mm in equivalent diameter
      • sandco_ : Mineral particles 0.5mm to 1.0mm in equivalent diameter
      • sandmed_: Mineral particles 0.25mm to 0.5mm in equivalent
      • sandfine_ : Mineral particles 0.10 to 0.25mm in equivalent diameter
      • sandvf_ : Mineral particles 0.05 to 0.10mm in equivalent diameter

3.2 Process

  • make buffer (400 m radius circles around well, which covers about 50 ha of land) for each well
source("./GitControlled/Codes/0_1_ls_packages.R")
source("./GitControlled/Codes/0_0_functions.R")
# /*===== point data for well location =====*/
data_w_LR_TB <- readRDS("./Shared/Data/WaterAnalysis/data_w_LR_TB.rds")

well_sf <- data_w_LR_TB %>%
  unique(.,by="wellid") %>%
  .[,.(wellid, longdd, latdd)] %>%
  # --- the geographic coordinate in the well data is NAD83 (espg=4269)--- #
  st_as_sf(., coords = c("longdd","latdd"), crs = 4269)


# /*===== make buffer around each well  =====*/
well_buffer_sf <- well_sf %>%
  #--- project to WGS UTM 14 (Cartesian 2D CS covering NE)---#
  st_transform(32614) %>%
  # --- make buffers--- #
  # 400 m radius circles around well, which covers about 50 ha of land
  st_buffer(., dist = 400)
 
ggplot()+
  geom_sf(data=well_sf, size=0.5) +
  geom_sf(data=well_buffer_sf, color="red", fill=NA)

well_buffer_sp <- as(well_buffer_sf, "Spatial")
  • Then, apply get_ssurgo_props() for each well buffer.
get_ssurgo_props <- function(field, vars, summarize = FALSE) {
  # field = as(filter(well_buffer_sf, wellid==112576), "Spatial")
  # field = as(filter(well_buffer_sf, wellid==112591), "Spatial")
  # field = as(filter(well_buffer_sf, wellid==231010), "Spatial")
  # field= as(well_buffer_sf[5,], "Spatial")
  # vars = soil_var
  # summarize = TRUE
  # Get SSURGO mukeys for polygon intersection
  ssurgo_geom <-
    SDA_spatialQuery(
      geom = field,
      what = "geom",
      db = "SSURGO",
      geomIntersection = TRUE
    ) %>%
    st_as_sf() %>%
    dplyr::mutate(
      # area = as.numeric(st_area(.)), #this causes an error for some polygons
      area = .$area_ac,
      area_weight = area / sum(area)
    )


  # --- check --- #
  # ggplot()+geom_sf(data=ssurgo_geom, aes(fill=factor(mukey)))
  # Looks like, in one buffer (surrounding a well), there are multiple polygon data where each of them are associated with a unique mukey 
  # `area_weight` indicates a portion of each polygon to entire butter

  # Get soil properties for each mukey
  mukeydata <-
    get_SDA_property(
      property = vars,
      method = "Weighted Average",
      mukeys = ssurgo_geom$mukey,
      top_depth = 0,
      bottom_depth = 150
    )

  ssurgo_data <- left_join(ssurgo_geom, mukeydata, by = "mukey")
  
  if (summarize == TRUE) {
    ssurgo_data_sum <-
      ssurgo_data %>%
      data.table() %>%
      .[,
        lapply(.SD, weighted.mean, w = area_weight, na.rm = TRUE),
        .SDcols = vars
      ]
    return(ssurgo_data_sum)
  } else {
    return(ssurgo_data)
  }
}

3.2.1 Example: SSURGO polygon data for an example field

library(soilDB)
library(aqp)
library(sp)

ex_well_bf <- well_buffer_sf[5,]
ex_well_bf_sp <- as(ex_well_bf, "Spatial")

ssurgo_geom <-
    SDA_spatialQuery(
      geom = ex_well_bf_sp,
      what = "geom",
      db = "SSURGO",
      geomIntersection = TRUE
    ) %>%
    st_as_sf()%>%
    dplyr::mutate(
      area = as.numeric(st_area(.)),
      area_weight = area / sum(area)
    )

# /*===== visualization =====*/
ggplot()+
  geom_sf(data=ssurgo_geom, aes(fill=factor(mukey)), color="blue", alpha=0.5)+
  geom_sf(data= well_sf[5,], fill=NA)

3.2.2 Download SSURGO data

  • This is the code I used to download SSURGO data for my data.
# === Target variables === #
soil_var <- c(
  # --- "Horizon" --- #
  "sandtotal_r", # a fraction of sand (0.05mm to 2.0mm)(%)
  "claytotal_r", # a fraction of clay (%)
  "silttotal_r", # a faction of silt (%) 
  "ksat_r", # hydraulic conductivity (um/m)
  "awc_r", # available water capacity (cm/cm)
  # --- "Component" --- #
  "slope_r" # The difference in elevation between two points, expressed as a percentage of the distance between those points. (SSM), "component"
  )

# each filed (buffer) data is passed on get_ssurgo_props()
res_ssurgo <- lapply(
  seq_len(nrow(well_buffer_sp)), 
  function(x) {
      print(paste0("working on wellid: ", well_buffer_sp[x, ]$wellid, " field"))
    # x=1
      get_ssurgo_props(
        field = well_buffer_sp[x,],
        vars = soil_var,
        summarize = TRUE
      ) %>%
      .[, wellid := well_buffer_sp[x, ]$wellid]
  }) %>%
  bind_rows()